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An infinite range spin glass like model for granular systems is introduced and studied through 
the replica mean field formahsm. Equilibrium, density dependent properties under vibration and 
gravity are obtained. 
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Handling of granular material is present in many agri- 
cultural and industrial processes and several fundamental 
practical problems are still unsolved. Besides that, their 
unusual static and flow properties 0] offer a challenging 
problem from the theoretical point of view, and despite 
the huge effort that has been devoted in recent years, 
are far from being fully understood. Since thermal en- 
ergy plays no role here, excitations can be achieved by 
externally shaking or shearing the system, enabling them 
to wander through the many microscopic configurations 
available for a fixed macroscopic density. Under vibra- 
tion, a multitude of fascinating phenomena show up, like 
heap formation, convection cells, size or shape segrega- 
tion, surface waves, etc (see ||l[ for a review and refer- 
ences). Another effect is the logarithmically slow rate of 
the density increase as the system suffers a sequence of 
taps This slow relaxation phenomenon under per- 
turbations, signalling complex cooperative movements of 
the particles, resembles the one found in systems with 
many metastable states like glasses and spin glasses. 

The analogy between glassy and granular behavior has 
been suggested some time ago and stressed recently 
along with the role of geometric frustration |^ . This frus- 
tration arises from the excluded volume of the grains, im- 
posing restrictions on their relative positions. For glass 
forming liquids, a simple frustrated lattice gas (FLG) has 
been introduced j|] that takes into account these steric 
effects and bridges complex fluids (glasses) and complex 
magnets (spin glasses). An inflnite range version [^||^] 
has also been studied in the framework of replica theory, 
yielding a very rich phase diagram. In order to apply this 
model for granular systems, Nicodemi et al ||] introduced 
the effects of gravity in the model, studying a tilted 2d 
lattice while applying a sequence of taps, the particles 
following a diffusion-like Monte Carlo dynamics. Among 
several interesting properties, analogous to those found 
in real experiments, the system displays an inverse loga- 
rithmic compactation behavior, reversible-irreversible cy- 
cles, aging, as well as a localization transition, signalled 
by a zero diffusion constant, in which the particles get 
trapped into dynamical local cages. This transition point 
seems to correspond to the Reynolds (or dilatancy) tran- 
sition observed in real systems. A lattice model has sev- 



eral advantages. From the computational point of view, 
the simulation cost is much lower than molecular dynam- 
ics. From the theoretical side, besides the amenability to 
theoretical investigation explored here, it helps in grasp- 
ing the fundamental concepts involved responsible for the 
complex observed phenomena. 

In the same spirit as Sherrington and Kirkpatrick |^ 
introduced an exactly solvable version of the Edwards- 
Anderson model pc| ] , here we introduce an infinite range 
version of the FLG considering L layers of N sites, con- 
nections being only between nearest neighbors layers (see 
figj^). Each site may be occupied by a particle [ui = 0, 1) 
having, for simplicity, only two possible spatial orienta- 
tion. Si ~ ±1. Although the grains can assume several 
spatial orientations (usually infinite), here we take the 
simplest case. The steric effects imposed by neighboring 
particles are felt as restrictions on the particle orienta- 
tion, what is included in the hamiltonian as quenched, 
gaussian distributed lattice bonds jf^ with jf^ = and 

{JfjY = J/N. Although the geometric frustration on 
the internal degrees of freedom should be considered as 
annealed at low densities and almost quenched at high 
ones, the quenched approximation is sufficiently good as 
can be seen from the 2d results Each layer has 

its own chemical potential satisfying /i^+i > (count- 
ing from top to bottom) that accounts for the effect of 
gravity and = gi/ L in order to have a constant force. 
Thus, we consider the following Hamiltonian: 
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The parameter K = —l-\-K' may tune the repul- 
sive/attractive {K' negative/positive) interaction be- 
tween particles H and may be important to treat wet 
powders. Here, as in the original model, we consider 
K = -J ^ -I {K' ^ Q). Notice that the value -1 ap- 
pears originally in order to recover, in the limit J 3> 1, 
the Frustrated Percolation (FP) Q constraint of only al- 
lowing non frustrated loops to be fully occupied. In this 
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limit, two given neighboring sites, only could be occu- 
pied {riinj = 1) if the corresponding orientations satis- 
fied the local disorder, JijSiSj = 1, otherwise at least 
one site should be empty {niTij =0). ft is also impor- 
tant to point out that by changing the value of K, sev- 
eral new qualitatively different phases do appear, in both 
frustrated or not versions of the Blumc-Emery-Griffiths 
(BEG) model (see and references therein). In the limit 
where all sites are occupied, we get a layered version of 
the SK model (L-SK). This model has a continuous tran- 
sition from a spin glass {qi ^ 0) phase to a paramagnetic 
one {q/, = 0,Vf) at ^ = l/^/2^c where Xc is the low- 
est positive root of a polynomial recursively obtained by 
P^^ix) = Pi^T^M - x'^Pi^^ix) with P(f ^ = Pf ^ = 1. 



-/3(pi£ + A^o) 
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A similar result |12| has been foimd for a modified mean- 
field version of the McCoy- Wu model This critical 
temperature approaches unity as L — s- cx). 
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The saddle point equations are rg = {qi-i + q£+i)/2 and 
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where go = = <lL+\ = d-L+i = 0. A global order 
parameter may be introduced as Q = L^^ qi. 

There is a critical temperature above which all qi 
are zero, that is, Q = 0. Figure ^ shows, for L = 100, 
the critical temperature as a function of density for sev- 
eral values of g. Notice that the bigger the density, the 
stronger should be the vibration in order to get a fluid 
state, similar to what happens for fixed density and in- 
creasing gravity. When g = (no gravity), de = p (V£) 
and Tc = pT^^ where T^^ is the critical temperature 
of the L-layered SK model (see below). On the other 
limit, when ^ — > oo (strong gravity), in analogy with 
the layered SK model, the critical temperature is related 
with the smallest positive root x* of a given polyno- 
mial, Tc — l/\/2a;*, depending on the density p. For 
< p < 1/L, we have Tc = because all particles occupy 
the lowest layer and do not interact. For 1/L < p < 2/L, 
some sites occupy the second lowest layer and the critical 
temperature is obtained with the root of P2{x). In gen- 
eral, for {n — l)/L < p < n/L, the relevant polynomial 
is Pn{x). These polynomials are obtained recursively by 



FIG. 1. System architecture with L layers, connections be- Pe{x) = 
ing only between nearest neighbors layers. 
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In evaluating the free energy, the trace is restricted to 
states with a given density p. Using standard techniques 
for dealing with disordered systems and assuming 
replica symmetry, the free energy reads 
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where Vz — (iz/V27rexp(— z^/2) and the temperature 
T = is a measure of the vibration imposed to the 
system. The order parameters are a diluted Edwards- 
Anderson = (^S°'^n°-^ S^^v!'^) and the density = 
(rf"^^ while /io accounts for the constraint p = L^^ dg 
and 



[qi+i + qe-i) 
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1 {di+i + df-i) 



with Pq{x) — Pi{x) = 1. In this limit, the critical tem- 
perature goes to unity for large L. 




FIG. 2. Phase diagram T versus p for L — 100 and g from 
(straight line) to 5 (top curve), showing the disordered phase 
[qi = 0,V/) and a spin glass phase {qt / 0). 
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The density and q profiles are shown in fig. g as a 
function of the system height £/L. Interestingly enough, 
the density profile shows, even at low temperature, that 
frustration effects are important in preventing a close 
packed configuration, signaled by a density lower than 
one. Its important to stress that, in analogy with the 
FLG 1^ , there are two regimes of densities depending on 
the gravity: one is the L-SK regime for large g with the 
particles settling in the lowest possible layers disregard 
the geometric effects and the second and most interesting 
one, shown here, where the steric effects become impor- 
tant. Although the g-profile seems to vanish above a 
given height, it actually does not do so for finite L and 
as we approach the continuous limit (L — > oo), a new 
transition settles down. From experiments and molec- 
ular dynamics simulation |16^, as a function of the 
vibration, the system passes from a solid-like behavior to 
a fluid-solid coexistence region. After this onset of flu- 
idization, the top layers become fluid-like, the particles 
having a great mobility. The interface between those re- 
gions decreases its height as the vibration increases. Al- 
though in this context a fluid regime means a situation 
where the particles present translational mobility, here 
we consider the case where the mobility is orientational, 
that is, a fluid layer will be one having qe — 0. How- 
ever, in some sense, both are a measure of the amount 
of geometrical constraints imposed on the particles by 
their neighbors. Moreover, when simulating the 3d sys- 
tem in absence of gravity [||, it can be seen that at the 
glass transition, the diffusion coefficient does vanish, cor- 
roborating our use of the term. For increasing values of 
L we can extrapolate the numerical results and find the 
critical layer, ic{T), separating the fluid region {qi = 0) 
and the solid one {qe ^ 0). The temperature where this 
first happens, signalling the onset of fluidization, is de- 
noted by Tp and the transition is discontinuous from Tp 
up to Ttrci while continuous for Ttrc < T < Tc- Notice 
that below Tp all non-empty layers have qe =^ 0. The 
point Ttrc is reminiscent of the tricritical point found in 
the disordered BEG model (remember that we have a 
varying chemical potential in the vertical axis). These 
information is summarized in figure ^. Notice that al- 
though in the figure ^ the top most layers have qe ~ 
below Tp, they are empty. Moreover, the density pro- 
file has been measured both in real experiments |T5| as 
well as in simulations, showing that in the steady state, 
independently of the phase of the up and down motion 
of the heap, the density profile is always preserved. This 
fact, indicating that configurational properties may be 
obtained through appropriate averages, independently of 
the dynamics, supports our equilibrium results and has 
also led to thermodynamic theories of granular systems 
p^ . Also, as the granular temperature increases, there 
is an elevation of the center of mass of the system [0 , 
given by hcM = {pL)~^ J2e^^i ^^'^ ^^^^ is ^ function of 
both temperature and gravity. 




FIG. 3. Density profile as a function of j/ = l/L for 
L — 1000, p — 0.5, g ~ 1 and several values of T. At high 
temperatures we recover the original mean field FLG as the 
density becomes uniform. Inset: profile of qt. Although the 
curve seems to go to zero above a certain layer, it actually 
does not do so for finite L. The onset of fluidization arrives 
at once for every layer. 




FIG. 4. Transition line 1 — £c/L showing the higher, flu- 
idized layers and the lower, solid {qe ^ 0) ones for g — 1, 
p — 0.5 and large L. The onset of fluidization occurs at 
Tp ~ 0.11 while the system is completely fluidized above 
Tc — 0.66. For T < Tp we have qe ^ for every non empty 
layer and since there is no longer a solid-fluid interface, just 
the material surface, we do not show the line. The transition 
is first order for Tp < T < Trc — 0.24 and continuous up to 



In conclusion, we introduced an infinite range version 
of a frustrated lattice gas model |^,|| for granular sys- 
tems and applied, to our knowledge for the first time, 
the replica formalism to these systems. In this mean 
field version, we are able to study stationary properties, 
obtaining the vibration, density and gravity dependent 
phase digram as well as information on the density pro- 
file and the onset of fluidization. 
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There is a multitude of possibilities that still remain to 
be explored. First of all, the stability of the replica sym- 
metric solution employed here and the behavior of the 
response functions, like the compressibility, for example. 
Also, the effects of segregation on polydisperse mixtures 
observed experimentally and in the 2d model , can 
be included by allowing different degrees of frustration 
for each type of particle in the Hamiltonian |19|. Results 
on the 2d model show that the out of equilibrium dy- 
namics, when the system is subject to small, continuous 
shaking [M, present aging in the two times correlation 
function C(t, t^) for the bulk density, remanescent of the 
glassy nature of the model. In view of this, it would be 
extremely interesting to study in detail the mean-field 
equations for the dynamics ||22| . Different coupling dis- 
tributions, e.g. bimodal, may be studied by means of the 
TAP formalism j2^, and in the limit of large L, the den- 
sity and q profile should obey a set of coupled differential 
equations, from which some analytical results may be ob- 
tained. From the simulational point of view, research is 
in course to study how compactation and segregation are 
affected by interpolating from the 2d model to the mean 
field case, along with the effects of attraction {K > — 1) 
between the particles, as in wet powders. Besides that, 
the fluidization transition may be detected on the simu- 
lation by measuring the mean squared displacement for 
each particle, Rlf{t) (instead of the system average) and 
comparing with its mean height. 
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